clear
t = 0:0.001*pi:2*pi;
z = exp(i*t);

%ABM_1
fz1 = z - 1;
plot(fz1, '-b', 'LineWidth', 2)
axis([-2.5, 0.5, -1.3, 1.7])
title('RAS plots for ABM with order 1,2,3')
hold on
num = 0;
pole = 0;
kp = -1;
r1 = roots([1, -1 - kp]);
k = 1;
for i= 1:k
    if norm(r1(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r1(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l1 = 'RAS of ABM with order 1, is bounded';
else
    l1 = 'RAS of ABM wit h order 1, is unbounded';
end

%ABM_2
fz2 = (z.*z - z)./(3/2*z - 1/2);
plot(fz2, '-r', 'LineWidth', 2)
num = 0;
pole = 0;
kp = -0.5;
r2 = roots([1, -1 - kp*3/2, kp*1/2]);
k = 2;
for i= 1:k
    if norm(r2(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r2(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l2 = 'RAS of ABM with order 2, is bounded';
else
    l2 = 'RAS of ABM with order 2, is unbounded';
end


%ABM_3
fz3 = (z.*z.*z - z.*z)./(5/12 - 16/12*z + 23/12*z.*z);
plot(fz3, '-y', 'LineWidth', 2)
num = 0;
pole = 0;
kp = -0.3;
r3 = roots([1, -1 - kp*23/12, kp*16/12, -kp*5/12]);
k = 1;
for i= 1:k
    if norm(r3(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r3(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l3 = 'RAS of ABM with order 3, is bounded';
else
    l3 = 'RAS of ABM with order 3, is unbounded';
end
legend(l1, l2, l3)
print(gcf, '-depsc','RAS_ABM_123.eps')
hold off

%ABM_4
fz4 = (z.*z.*z.*(z - 1))./(-9/24 + 37/24*z - 59/24*z.*z + 55/24*z.*z.*z);
x41 = [0];
y41 = [0];
x42 = [0];
y42 = [0];
x43 = [0];
y43 = [0];
for i = 1:2001
    if imag(fz4(i))< 0.4311 && imag(fz4(i)) > -0.4311
       x41 = [x41 real(fz4(i))];
       y41 = [y41 imag(fz4(i))];
    elseif imag(fz4(i) > 0.4311)
       x42 = [x42 real(fz4(i))];
       y42 = [y42 imag(fz4(i))];
    else 
       x43 = [x43 real(fz4(i))];
       y43 = [y43 imag(fz4(i))];
    end
end
plot(fz4, '-b')
title('RAS plot of ABM with order 4')
axis([-1, 1, -1, 1])
hold on
kp = -0.1;
r4 = roots([1, -1-kp*55/24, kp*59/24, -kp*37/24, kp*9/24]);
k = 4;
num = 0;
pole = 0;
for i= 1:k
    if norm(r4(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r4(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
     fill(x41, y41, '-y')
end
num = 0;
pole = 0;
kp = 0.1 + 0.6i;
r4 = roots([1, -1-kp*55/24, kp*59/24, -kp*37/24, kp*9/24]);
for i= 1:k
    if norm(r4(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r4(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
     fill(x42, y42, '-y')
end
num = 0;
pole = 0;
kp = 0.1 - 0.6i;
r4 = roots([1, -1-kp*55/24, kp*59/24, -kp*37/24, kp*9/24]);
for i= 1:k
    if norm(r4(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r4(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
     fill(x43, y43, '-y')
end
print(gcf, '-depsc','RAS_ABM_4.eps')  
hold off

clear
t = 0:0.001*pi:2*pi;
z = exp(i*t);
%AMM_3
fz3 = z.*(z-1)./(5/12*z.*z + 8/12*z - 1/12);
plot(fz3, '-b', 'LineWidth', 2)
axis([-7, 1, -3.2, 5.2])
title('RAS plots for AMM with order 3,4,5')
hold on
num = 0;
pole = 0;
kp = -3;
r3 = roots([1 - kp*5/12, -1 - kp*8/12, kp*1/12]);
k = 2;
for i= 1:k
    if norm(r3(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r3(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l3 = 'RAS of AMM with order 3, is bounded';
else
    l3 = 'RAS of AMM wit h order 3, is unbounded';
end

%AMM_4
fz4 = z.*z.*(z - 1)./(9/24*z.*z.*z + 19/24*z.*z - 5/24*z + 1/24);
plot(fz4, '-r', 'LineWidth', 2)
hold on
num = 0;
pole = 0;
kp = -1;
r4 = roots([1 - kp*9/24, -1 - kp*19/24, kp*5/24, -kp*1/24]);
k = 3;
for i= 1:k
    if norm(r4(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r4(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l4 = 'RAS of AMM with order 4, is bounded';
else
    l4 = 'RAS of AMM wit h order 4, is unbounded';
end

%ABM_5
fz5 = z.*z.*z.*(z - 1)./(251/720*z.*z.*z.*z + 646/720*z.*z.*z - 264/720*z.*z + 106/720*z - 19/720);
plot(fz5, '-y', 'LineWidth', 2)
hold on
num = 0;
pole = 0;
kp = -1;
r5 = roots([1 - kp*251/720, -1 - kp*646/720, kp*264/720, -kp*106/720, kp*19/720]);
k = 4;
for i= 1:k
    if norm(r5(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r5(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l5 = 'RAS of AMM with order 5, is bounded';
else
    l5 = 'RAS of AMM wit h order 5, is unbounded';
end
legend(l3, l4, l5)
print(gcf, '-depsc','RAS_AMM_345.eps')
hold off

clear
t = 0:0.001*pi:2*pi;
z = exp(i*t);
%BDF_1
fz1 = (z-1)./(z);
plot(fz1, '-b', 'LineWidth', 2)
axis([-2, 12, -7, 10.5])
title('RAS plots for BDF with order 1,2,3,4')
hold on
num = 0;
pole = 0;
kp = 0.5;
r1 = roots([1 - kp*1, -1]);
k = 1;
for i= 1:k
    if norm(r1(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r1(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l1 = 'RAS of BDF with order 1, is bounded';
else
    l1 = 'RAS of BDF with order 1, is unbounded';
end

%BDF_2
fz2 = (z.*z - 4/3*z + 1/3*z)./(2/3*z.*z);
plot(fz2, '-r', 'LineWidth', 2)
hold on
num = 0;
pole = 0;
kp = 2;
r2 = roots([1 - kp*2/3, -4/3, 1/3]);
k = 2;
for i= 1:k
    if norm(r2(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r2(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l2 = 'RAS of BDF with order 2, is bounded';
else
    l2 = 'RAS of BDF with order 2, is unbounded';
end

%ABM_3
fz3 = (z.*z.*z - 18/11*z.*z + 9/11*z - 2/11)./(6/11*z.*z.*z);
plot(fz3, '-y', 'LineWidth', 2)
hold on
num = 0;
pole = 0;
kp = 4;
r3 = roots([1 - kp*6/11, -18/11, 9/11, -2/11]);
k = 3;
for i= 1:k
    if norm(r3(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r3(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l3 = 'RAS of BDF with order 3, is bounded';
else
    l3 = 'RAS of BDF with order 3, is unbounded';
end

%ABM_4
fz4 = (z.*z.*z.*z - 48/25*z.*z.*z + 36/25*z.*z - 16/25*z + 3/25)./(12/25*z.*z.*z.*z);
plot(fz4, '-m', 'LineWidth', 2)
hold on
num = 0;
pole = 0;
kp = 4;
r4 = roots([1 - kp*12/25, -48/25, 36/25, -16/25, 3/25]);
k = 4;
for i= 1:k
    if norm(r4(i), 2) == 1
        num = num + 1;
        if(num > 1)
            pole = 1;
        end
    end
    if norm(r4(i), 2) > 1
        pole = 1;
    end
end
if pole == 0
    l4 = 'RAS of BDF with order 4, is bounded';
else
    l4 = 'RAS of BDF with order 4, is unbounded';
end
legend(l1, l2, l3, l4)
print(gcf, '-depsc','RAS_BDF_1234.eps')
hold off